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1. Introduction 

It is well known and appreciated that the accretion of ice on lifting surfaces can significantly degrade the 
performance and handling characteristics of affected aircraft [1], Less well understood, however, are the de- 
tails of the complex flow fields associated with these configurations. These flow fields may exhibit extensive 
regions of unsteady separated flow depending on the ice shape. Typically, ice accretions on straight wings 
arc categorized as “rime,” “horn,” or “spanwise ridge,” which arc more or less two-dimensional in nature, 
or “roughness,” which is the most three-dimensional of the shapes [1], The flow field for a horn ice accre- 
tion, which is the subject of this effort, is dominated by a separation bubble downstream of the horn. The 
flow separates because it is unable to recover the necessary pressure in the boundary layer after it expands 
around the horn. The ensuing free shear layer rolls up to form vortices that arc convected downstream. The 
development of this unsteady flow field is critical to the resulting flow near maximum lift. The separation 
bubble has a large global effect and is similar to the long bubble defined by Tani [2], 

Many researchers have turned to computational fluid dynamics (CFD) simulations to investigate iced 
wing flow fields. Beginning with Potapczuk [3], numerous results from icing effects studies using the 
Reynolds averaged Navier-Stokes (RANS) equations have been reported in the literature. Chung, et al. [4], 
analyzed the flow around the ice contaminated wing surfaces of a turbo-prop aircraft. Their two- and three- 
dimensional analyses of a wing with a spanwise ridge accretion employed a steady RANS approach and 
were performed to provide insight into the aerodynamics that may have led to a loss of control of the 
aircraft. Reported discrepancies between the two- and three-dimensional results were attributed to a lack 
of grid resolution for the three-dimensional cases considered. Dunn and Loth [5] studied the effects of 
simulated-spanwise-ice shapes on the aerodynamic performance of two-dimensional wing sections using 
a RANS solver developed for unstructured meshes. A forward-facing quarter round was employed as the 
simulated ice shape. Additionally, they employed a solution adaptive mesh in which the mesh was refined 
in the separated shear layer. Their results showed that, even with adaptive mesh refinement, the pressure 
recovery was not correctly predicted. This discrepancy was attributed to the inability of the steady RANS 
simulation to model the experimentally observed unsteadiness that occurs in flow fields of this type. Pan, 
Loth, and Bragg [6] performed RANS simulations for flow about airfoils with ridge ice shapes and leading 
edge ice shapes and compared the results to experimental data. Favorable comparisons between predicted 
force data and experiments were obtained up to, but not including, the stall condition, which is dominated 
by unsteady flow. 

Recently, researchers have begun employing unsteady simulations to investigate separated flow fields. 
Candidate approaches for computing such unsteady flow fields include unsteady RANS (URANS), direct 
numerical simulation (DNS), large-eddy simulation, (LES), and detached eddy simulation (DES). Spalart [7] 
contends that URANS is both “ambiguous and flawed” and that “its quantitative performance can be quite 
poor.” DNS, while attractive from the standpoint that no modeling is required, and LES are currently im- 
practical for realistic configurations due to their computational expense. DES, however, appeal's to have 
the potential to represent a viable near-term approach. DES is a hybrid RANS/LES approach in which a 
RANS turbulence model is employed in attached thin shear layers and an LES-type model is employed in 
regions away from the wall. This approach exploits the ability of RANS simulations to efficiently model 
high Reynolds number attached boundary layer flows and the ability of LES to model geometry dependent, 
unsteady three-dimensional flows [8]. DES has been employed in numerous simulations of complex flow 
fields [9-1 1]. DES has also been applied to the icing effects problem. Kumar and Loth [12] employed a DES 
technique to predict the unsteady flow around an NLF0414 airfoil section with a synthetic ice shape located 
on the upper surface at 3.4% chord. They also reported on computations for a rectangular three-dimensional 
wing using the same section. Their three-dimensional results showed some improvement relative to the 
steady RANS computations. Pan and Loth [13] performed a DES for an NACA23012 airfoil with forward- 
facing quarter-round simulated ice accretion. Again, some improvement between the predicted results and 
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experimental data was reported. 

We previously reported results from steady RANS simulations for a horn ice shape/wing configura- 
tion [14, 15]. The basic ice shape considered was a 22.5-minute glaze ice accretion on a GLC305 airfoil, 
which is denoted as the two-dimensional 944-ice shape [16]. The airfoil section/ice shape was extruded to 
form a rectangular planform wing with an aspect ratio of unity. The flow fields associated with this con- 
figuration have been studied experimentally and the results reported in [17, 18]. These data include force 
measurements and pressure measurements as well as mean velocity measurements and RMS fluctuations of 
the streamwise and transverse velocity components obtained using a split-hot-film probe. The model was 
mounted between the tunnel walls, so no tip effects were considered. We compared steady RANS results 
computed using the noncommercial version of the Cobalt^ code [19] (now called AVUS) and the Spalart- 
Allmaras turbulence model [20] with the experimental data. In general, the steady RANS computations 
underpredicted both the lift and drag as the angle of attack was increased. Additionally, the extent of the 
separated region was overpredicted even at lower angles of attack. As noted above for other RANS compu- 
tations, the pressure recovery was not in agreement with the experimental data. In a related effort, Chi, et 
al. [21] reported on the results of a series of two-dimensional steady RANS simulations for the GLC305/944- 
ice shape, as well as a rime ice shape, using a variety of turbulence models. The two-dimensional results for 
the 944-ice shape were qualitatively very similar to those reported in Ref. [14, 15] regardless of the turbu- 
lence model employed. Additionally, relatively fine meshes were employed to discretize the computational 
domain. Therefore, we believe the discrepancies reported in Ref. [21] were not related to mesh effects. 

In this report, we focus on activities associated with the first year of an effort to investigate the capability 
of existing flow simulation methodologies to predict the separated, highly unsteady flow fields associated 
with wings with significant horn-ice accretions. Here, we present the results of steady RANS computations 
and unsteady DES computations performed for the rectangular planform, extruded GLC305/944-ice shape 
wing. We first provide a brief discussion of the DES approach. We then provide a more concrete definition 
of the problem. Numerical considerations regarding how the mesh spacing and time step size are selected 
along with the mesh generation and solution procedures are discussed next. We present sample results 
from two-dimensional and three-dimensional computations. We then provide an assessment of the tools we 
employed and make a few recommendations regarding the roles of RANS and DES computations for icing 
effects studies. We conclude with a brief summary of the effort. 
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2. Detached Eddy Simulation 

DES was originally formulated as a modification to the Spalart-Allmaras turbulence model [8,20]. The 
basic idea is to employ a RANS-type turbulence model for thin, attached shear layers and an LES model 
in separated regions. This approach exploits the ability of RANS simulations to efficiently model attached 
boundary layer flows and the ability of LES to model unsteady, geometry-dependent, three-dimensional 
flow fields. RANS-type models perform well for attached flows and less well for separated flows, while 
LES models perform well for separated flows but arc prohibitively expensive to employ for computations 
involving thin boundary layers. We now provide a brief overview of the DES model as implemented in the 
Spalart-Allmaras turbulence model. More detailed discussions can be found in References [8, 10, 1 1, 22]. 

In the Spalart-Allmaras model, the wall destruction term is taken to be proportional to (D/d) 2 where 
v is a working variable related to the turbulent viscosity and d is the distance to the nearest wall. If local 
equilibrium occurs and the production and destruction terms are in balance, the eddy viscosity becomes 
proportional to Sd 2 where S is the local strain rate. In the Smagorinsky LES model [23], the subgrid-scale 
turbulent viscosity is proportional to the local strain rate and the local mesh spacing squared, i.e., .S' A 2 
where A = max (Ax, Ay, Az). Thus, if the distance to the wall d is replaced by the local mesh spacing 
A, the Spalart-Allmaras model will locally behave like a Smagorinsky LES model. To retain the RANS- 
type behavior in attached boundary layers, d is replaced by a new variable d = min (d, Cdes^) where 
Cdes is a constant. Note that when d « CdesA , the model is in a RANS mode and models the average 
properties of attached flow turbulence. When d » CdesA, the model is in LES mode and resolves 
eddies larger than some wavelength depending on the characteristics of the flow solved and the local mesh 
spacing. For unstructured meshes, the spacing A is taken to be the longest distance between the cell center 
and neighboring cell centers [22]. The implementation of the DES model has two immediate implications 
regarding the mesh: 

1. The mesh in the boundary layer should be highly anisotropic. In particular, the spacing along the 
surface should be larger than the local boundary layer thickness. This ensures that the model operates 
in a RANS mode near the boundary. 

2. There is no advantage to having an anisotropic mesh in the interior of the domain. The premise 
of LES is to filter out only eddies that arc statistically isotropic [7] so that equal resolution in all 
directions is reasonable. 

Other implications regarding the discretization of the domain are discussed in the Section 3.3. 
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3. Numerical Considerations 

The prediction of the aerodynamic characteristics of an iced wing is a complex problem that involves several 
steps. Given a geometric definition of a wing with an ice accretion, these steps must be followed: 

• Develop a representation of the surface that is suitable for generating a surface mesh. 

• Specify the artificial boundaries needed to define the computational domain, e.g., the outer boundary 
and the side boundaries. 

• Generate a mesh and specify boundary conditions on the bounding surfaces of the computational 
domain. 

• Generate a mesh in the interior of the computational domain. 

• Generate a flow solution. 

• Analyze the results. 

The following sections describe how above steps were accomplished in this effort. 

3. 1 Problem statement 

The specific airfoil section considered during this effort was the GLC305 airfoil section with the 22.5- 
minute glaze ice accretion which is denoted as the two-dimensional 944-ice shape [16]. The ice shape was 
extruded to form a rectangular planfomr wing with an aspect ratio of unity to match the configuration used 
in the test program [17, 18] and is referred to here as the “extruded wing.” Since the wing employed for 
the test program was mounted between walls, the geometry modeled in this study did not include wing 
tips. Numerical solutions were compared with experimental data for the following conditions: M=0.12, 
Re/L=3.8xl0 6 /m which, with a chord length of 0.9144m, yields Re=3.5xl0 6 , and angles of attack of 0°, 
2°, 4° and 6°. The extruded 944-ice shape cases considered here correspond to Run 41 in the experimental 
data [17, 18], 

3.2 Geometry modeling 

The GLC305/944-ice shape definition was used as input to ICEG2D [24] , a two-dimensional tool that 
automates geometry modeling and mesh generation for ice accretion predictions. ICEG2D redistributes 
points on the upper and lower surfaces of the defined airfoil using a curvature-based equidistribution al- 
gorithm to ensure a sufficient number of points arc employed in regions of the airfoil surface with high 
curvature. A structured surface mesh was then extruded from the section definition. This surface mesh was 
then converted to a meshing -ready NURBS representation using the mesh generation software GUM-B [25]. 

3.3 Mesh Generation: Estimated DES zones for iced wing 

It is well known that the quality of the CFD solution depends on an appropriate mesh with sufficient 
mesh density in regions of high flow gradients. In an unsteady flow dominated by convecting vortices, these 
vortices must be resolved in addition to attached boundary layers. Indeed, as noted by S pal art [7], “ DES 
compounds the gridding difficulty by incorporating both types of turbulence treatment in the same field.” 
Spalart provides guidance as to how a domain should be discretized for DES modeling. He identifies several 
different regions with different meshing requirements: 

• An Euler region (ER) is a region that is free of turbulence and vorticity unless it is penetrated by a 
shockwave. The ER covers most of the domain and an isotropic mesh can be employed in this region. 

• A RANS region (RR) is primarily composed of the boundary layer where there is no LES content. 
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Figure 1: Estimated DES zones for GLC305 airfoil section with two-dimensional 944-ice shape 
(extruded wing) 

- The viscous region (VR) is located within the RANS region and its meshing requirements are 
similar to those of standard RANS computations. 

- The RANS outer region (OR) should be discretized with a mesh in which the mesh spacing 
in the direction normal to the wall should not exceed once-tenth the thickness of the boundary 
layer. This is primarily an issue associated with numerical robustness. In practice, this rule is 
often violated. 

• An LES region (LR) contains vorticity and turbulence but is not a boundary layer. Significant LES 
content is present in an LR. 

- The focus region (FR) is the region near the body in which the separated turbulence must be 
well resolved. The mesh should be isotropic in the FR since the LES mode filters out eddies 
that are statistically isotropic. 

- The departure region (DR) is a transitional region between the FR and the ER. The DR does 
not need to be resolved as well as the FR. 

Here, regions are not distinguished by different equations being applied but by different priorities in the 
mesh spacing. An efficient mesh for any external flow can be designed with these concepts in mind, but not 
all are strict requirements. Figure 1 depicts the estimated DES zones discussed above as applied to the flow 
surrounding the 944-ice shape. 

3.4 Selection of mesh spacing in focus region 

The procedure described below follows that suggested by Spalart [7]. According to S pal art. a well- 
adjusted subgrid-scale model should allow energy cascade to the smallest eddies that can be tracked on the 
mesh. Therefore, for most CFD solvers, an eddy with a wavelength of A = 5Ao, where A 0 is the local 
mesh spacing, will be active even though it cannot be highly accurate because it lacks the energy cascade to 
smaller eddies, and is under the influence of eddy viscosity instead. In other words, if we want to adequately 
resolve an eddy with wavelength A, the mesh spacing should be Ao = A/ 5. In the baseline mesh employed 
for DES computations, an eddy with a wavelength of 5% of the chord was selected, resulting in mesh 
spacing in the focus region of 1% of the chord. This wavelength was chosen because the height of the horn 
is approximately 5% of the chord. Therefore, vortices with approximately the same spatial extent as the 
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horn height can be resolved in the separated region downstream of the horn. As paid of the comparison of 
the DES results with experimental data, the DES mesh was refined by a factor of two in the focus region. 

3.5 Selection of time-step size 

Selecting the time-step size for DES computations is challenging. Spalart [20] suggests employing a 
local CFL number (based on the local flow velocity, the local mesh spacing, and the time step size) of unity 
in the focus regions (FR), that is where Uo is the maximum flow velocity in the region. Computed 

RANS results indicate that the maximum flow velocity over the domain is approximately 40-50% greater 
than the freestream velocity. For the problem considered here with a freestream velocity of approximately 
41.1m/s and a wing chord of 0.9144m, the computation for the time step size yields At=0. 1 5ms for the 
baseline DES mesh. For the refined mesh, zD=0. 075ms. However, as Spalart notes [7], steps a factor of 
3/2 or even 2 away in either direction from this estimate cannot be considered as “incorrect.” Unfortunately, 
tests with different time steps rarely give any strong indications toward an optimal value [7]. Thus, some 
ambiguity remains in the selection of the time step size. Solutions were therefore obtained on the baseline 
DES mesh using the smaller time step to test the effects of reducing the time step size for a fixed mesh. 

3.6 Unstructured mesh generation using Sol id Mesh and VGridns/GridTool 

Three different meshes were generated for the extruded wing configuration: a relatively coarse mesh 
that was used solely for RANS computations, a baseline DES mesh, and a refined DES mesh. The coarse 
RANS mesh was generated using SolidMesh. SolidMesh is an interface to the unstructured surface and 
volume mesh generation software AFLR2 and AFLR3 [26]. AFLR3 uses an advancing front algorithm to 
insert a point in the mesh. The point insertion is followed by a local reconnection to improve mesh quality. 
The baseline DES and refined DES meshes were generated using GridTool [27] and VGridns [28]. VGridns 
uses an advancing layer algorithm to generate a tetrahedral volume mesh. The near-body elements of these 
all-tetrahedral meshes were converted into prisms using the Blacksmith utility [29] thereby producing a 
mixed element hybrid mesh. SolidMesh produces a hybrid mesh with prisms and tetrahedra automatically. 
These hybrid meshes were employed because of their potential for improved efficiency and accuracy in 
comparison to unstructured tetrahedral meshes. The GridTool/VGridns combination was chosen for the 
DES meshes because it gives the user reasonably good control of the point spacing on the wing surface and 
the flow field through the use of “sources.” 

Figure 2 shows the source arrangement for the GLC305 wing with the 944-ice shape. Mesh refinement 
can be obtained simply by reducing the magnitude of the sources in appropriate regions. In this case, the 
mesh was locally refined by reducing the source strengths in the upper surface FR by a factor of two while 
holding the other source strengths constant. This was found to be a reasonably efficient mechanism for 
mesh refinement for this problem. Since the focus of this effort was the flow on the upper surface, the 
lower surface mesh was made relatively coarse in comparison. This is indicated by the relative sizes of the 
midchord source elements for the upper and lower surfaces in Figure 2. 

Figures 3(a), 3(b), and 3(c) show cross sections through the coarse RANS, baseline DES, and refined 
DES meshes respectively. These figures show cutting planes so that the line segments represent intersections 
of cell faces with the cutting plane. The connectivity of the mesh is evident and the tetrahedra/prism layers 
near the surface are clearly visible. Figure 3(a) shows the faceted surface represented by the triangular 
surface mesh on the extruded wing. Additionally, due to the manner in which the spacing is controlled in 
SolidMesh, there is a spanwise variation in the spacing. The volume mesh in the region downstream of 
the upper horn is relatively coarse. Refining the surface mesh reduces the faceting and improves the mesh 
density in the region where flow separation is anticipated, as shown in Figures 3(b) and 3(c). Also note 
the finer resolution in the focus region downstream of the horn for the refined mesh. These figures also 
illustrate the connection between the surface mesh and volume mesh. Because of the manner in which the 
mesh is generated, i.e., anisotropic tetrahedra/prisms transitioning to isotropic tetrahedra, the surface mesh 
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Figure 2: “Sources” used by GridTool/VGridns to generate mesh 


characteristics are propagated into the volume mesh. Therefore, in order to refine the mesh spacing in the 
domain, it is necessary to refine the mesh on the surface. 

The meshes shown in Figures 4(a) and 4(b) were used for the two-dimensional computations. Steady 
state two-dimensional RANS and unsteady two-dimensional DES solutions were obtained on the baseline 
mesh shown in Figure 4(a) and the refined mesh shown in Figure 4(b). The meshes employed for the two- 
dimensional simulations were generated by extracting the surface mesh in the symmetry plane from the 
three-dimensional meshes. 

Table 1 shows statistics for the five meshes employed to generate the results reported here. Notice that 
the localized mesh refinement in the refined DES mesh produces more than a 25% increase in the number 
of cells when compared to the baseline DES mesh. In some sense, this mesh refinement can be thought of 
as a manual, albeit crude, attempt at solution adaptive meshing. In all cases, the distance to the first point 
off the wall was defined so that an average y + value less than 0.5 was obtained. This value is well within 
the recommended values for the turbulence model employed. 

3.7 Flow solution 

The flow solver employed in this effort is the non-commercial version of Cobalteo [30], which is now 
called AVUS to reduce confusion with the commercial flow solver. The AVUS flow solver was designed for 
general unstructured meshes. It employs a nonlinear Riemann solver for the inviscid flux computations and 
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(a) Two-dimensional baseline mesh 



(b) Two-dimensional refined mesh 


Figure 4: Two-dimensional mesh for GLC305/944-ice shape 

can be run either in explicit or implicit mode. Second-order spatial accuracy is obtained using a linear least- 
squares reconstruction of the data. The detailed computational methodology used in AVUS is described in 
reference [19]. 

All computations were performed using second-order spatial accuracy. The RANS solutions reported 
here were obtained using first-order, implicit local time stepping and do not in any way represent time- 
accurate solutions. For the DES computations, the second-order temporal integration scheme was employed. 
It should be noted that a CFL number of unity in the focus region yields a CFL number on the order of 
100,000 in the boundary layer. Therefore, for stability reasons, it was necessary to use the fully-implicit 
integration scheme { 6 = 1 ) for the DES computations. Based on recommendations in the Cobalt^ user’s 
manual [29], two Newton iterations, each consisting of 30 Gauss-Seidel iterations, were employed per time 
step. Several turbulence models arc available including the Spalart-Allmaras one-equation model [20], 
which was used for the computations reported here. No transition point was specified, and the flow was 
assumed to be fully turbulent. A slip boundary condition was applied on the artificial side boundaries to 
reduce computational expense. 
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Table 1: Statistics for the meshes employed in this effort 



# Nodes 

# Faces 

# Cells 

Extruded Wing (coarse RANS mesh) 

890,958 

5,299,738 

2,258,447 

Extruded Wing (baseline DES mesh) 

1,828,711 

17,327,367 

8,341,019 

Extruded Wing (refined DES mesh) 

2,394,393 

22,088,508 

10,579,834 

Two-dimensional (baseline DES mesh) 

26,456 

48,775 

N/A 

Two-dimensional (refined DES mesh) 

28,656 

52,128 

N/A 
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4. Results 

This section presents detailed comparisons with experimental data [17, 18] of the steady RANS and unsteady 
DES results obtained for the extruded wing (GLC305/944-ice shape). In section 4.1, details of the experi- 
ments employed for validation arc presented. Comparisons of results obtained from the three-dimensional 
RANS simulations with experimental data for the extruded wing arc given in section 4.2. The results ob- 
tained from two-dimensional and three-dimensional DES results arc compared with experimental data in 
section 4.3. 

The RANS simulations were performed using an impulsive start and were continued until the normal 
force reached an approximate steady state, which typically required 15,000-16,000 global iterations. The 
RANS solutions were then used as initial conditions for the DES computations. The sudden change in the 
eddy viscosity (due to the reduction in the length scale when the DES computation is initiated) represents a 
non-physical transient that must be eliminated. Therefore, in all cases, DES computations were performed 
for 0.75s (a nondimensional time of 33.7), which corresponds to 5000 steps for At=0. 1 5ms and 10000 steps 
for At =0.07 5 ms to eliminate this transient. Data was then collected during an additional 0.75s for a total 
nondimensional time of 67.4 beyond the steady RANS solution. 

All cases reported in this effort were run on 64 processors on the EMPIRE cluster or the MAVERICK 
cluster at the ERC at Mississippi State University. The EMPIRE cluster is a supercomputer class cluster 
of workstations consisting of 1038 one GHz or better Pentium III processors each with one or more GB of 
RAM. The MAVERICK cluster is also a supercomputer class cluster comprised of 192 IBM x335 nodes. 
Each node contains dual 3.06 GHz Xeon processors and 2.5 GB of RAM. We now provide timing data for 
the runs performed using the EMPIRE cluster. Simulations performed on the MAVERICK cluster required 
approximately one-half of the time of those performed on the EMPIRE cluster. 

For the three-dimensional RANS solutions on the “baseline DES mesh,” 1 hour and 39 minutes were 
required for each 100 iterations on the EMPIRE cluster, whereas 2 hours and 23 minutes were required for 
each 100 iterations on the “refined DES mesh.” Therefore, a steady RANS computation for a single angle 
of attack required approximately 250 hours on the “baseline DES mesh” and 350 hours on the “refined DES 
mesh.” Here, hours refers to CPU hours. 

For the three-dimensional DES computations, 2 hours and 52 minutes were required for 100 time steps 
on the “baseline DES mesh.” Three hours and 52 minutes were required for each 100 time steps on the “re- 
fined DES mesh.” DES computations using the standard time step for the refined mesh, i.e., xU=0. 075ms, 
were performed for 20000 time steps and required a total of approximately 575 hours for the “baseline DES 
mesh” and 775 hours for the “refined DES mesh.” Thus, including the time required to compute the con- 
verged RANS simulations from which the DES computations were initiated, approximately 820 hours were 
required to compute the DES results for one angle of attack on the “baseline DES mesh.” Approximately 
1 125 hours were required to compute the DES results for a single angle of attack on the “refined DES mesh.” 

4. 1 Experimental details 

Addy et al. [17] performed icing effects studies for a 36-inch chord, two-dimensional business jet airfoil 
(GLC305) by conducting wind tunnel tests in NASA Langley Research Center’s LTPT. The GLC305 airfoil 
is designed for low transonic drag and has a maximum thickness to chord ratio of 8.7%. Four different types 
of ice shapes were considered for this study. Ice shapes were accreted on the GLC305 airfoil in the IRT 
over a range of icing conditions selected from the FAA FAR Paid 25-Appendix C [31]. The 22.5-minute 
glaze 944-ice shape was one of the ice shapes considered. The Reynolds numbers ranged from 3.0xl0 6 
to 10.5xl0 6 at a fixed Mach number of 0.12. The effects of Mach number variation (0.12 and 0.28) were 
investigated at constant Reynolds numbers of 6.0xl0 6 and 10.5xl0 6 . 

The airfoil model was supported horizontally across the width of the test section between two 40-inch 
diameter circular end-plates. These end-plates were flush mounted with the sidewalls and rotated to provide 
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the angle of attack adjustment. Each end-plate was equipped with a porous section for sidewall boundary 
layer control. The lift and pitching moment data were obtained from the integration of surface static pres- 
sures, while the drag coefficients were calculated using the standard momentum deficit method based on the 
pressures measured using a wake probe. Corrections to the integrated performance coefficients accounting 
for solid and wake blockage and streamline curvature were applied to the data during post processing. 

Figures 5(a) and 5(b) depict comparisons of the lift and drag coefficients, respectively, for the clean and 
the extruded wings. Figure 6 shows the pressure coefficient for the clean and extruded wings for an angle 
of attack range of 0° to 6° in 2° increments. The data shown here was obtained at a Mach number of 0.12 
and a Reynolds number of 3.5xl0 6 for the extruded wing case while the clean wing data was obtained at the 
same Mach number and a Reynolds number of 3.0xl0 6 . 

Broeren et al. [18] also carried out flow field measurements on the upper surface of the same model. 
Data were obtained at three different angles of attack preceding stall at Reynolds numbers of 3.5xl0 6 and 
6.0xl0 6 and Mach numbers of 0.12 and 0.21. Split-hot-film anemometry was used to measure the time- 
averaged flow velocities and its RMS fluctuations. Results confirmed the presence of a large separation 
bubble downstream of the ice shapes. Figure 7 shows experimental time-averaged stream wise velocity plots 
for the extruded wing at angles of attack of 0°, 4° and 6°. In general, these time-averaged velocity contour 
images provide a good overall illustration of the separated flow past the ice shape. From these images, we 
can see how the boundary layer separates near the tip of the ice horn. A significant reversed flow region is 
also formed below the separated shear layer. 
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(a) Experimental lift coefficients 



(b) Experimental drag coefficients 


Figure 5: Comparison of experimental lift and drag coefficients for extruded GLC305 airfoil 
section: with and without 944-ice accretion [17] 
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(c) 4° angle of attack 


(d) 6° angle of attack 


Figure 6: Comparison of experimental pressure coefficients for extruded GLC305 airfoil section: 
with and without 944-ice accretion [17] 
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Figure 7: Experimental u-velocity contours for the extruded wing at angles of attack of 0°, 4° and 
6° [18] 


NAS A/CR— 2004-2 13379 


15 


4.2 RANS results 

In this section, steady RANS results obtained for the extruded wing are compared with experimental 
data [17, 18]. RANS computations were performed on the meshes denoted as “coarse RANS,” “baseline 
DES,” and “refined DES,” shown in Figures 3(a), 3(b) and 3(c), respectively, at angles of attack of 0° to 6° 
in 2° increments. For angles of attack of 0° and 2°, 14,000 iterations were performed on the baseline DES 
mesh and refined DES mesh. For angles of attack of 4° and 6°, 16,000 iterations were performed on both 
the baseline and refined DES meshes. Results arc also included for the two-dimensional simulations for 
reference. 

Figure 8 shows the convergence histories for the extruded wing computation for an angle of attack of 6° 
using the baseline DES mesh. The lift force and axial force reached asymptotic values indicating a converged 
solution as shown in Figures 8(a) and 8(b), respectively. The 6° convergence history is representative of the 
computations for other angles of attack using this mesh. 

Figure 9 shows comparisons of the wing lift and drag coefficients computed from the steady RANS 
simulations with experimental data for an angle of attack range of 0° to 6° in 2° increments. As shown in 
Figure 9(a), the lift coefficient at 0° shows good agreement with experimental data for all cases, although 
the two-dimensional cases arc somewhat underpredicted. However, as the angle of attack increased, the 
deviation between the predicted lift and experimental values increased. The lift is underpredicted relative 
to the experimental data. Additionally, the premature break in the predicted lift curve slope, characteristic 
of the “near stall” behavior of the wing, suggests that the predictions arc overestimating the severity of 
the flow separation. As shown in Figure 9(b), the predicted drag coefficient shows reasonable agreement 
with the experimental data. This relatively good agreement (for drag) can be attributed to the fact that the 
drag is primarily composed of form or pressure drag for this configuration, i.e., the ice shape is a bluff 
body. However, the good agreement between the results predicted using the coarse RANS mesh and the 
experimental data appeal's to be somewhat fortuitous. The agreement is not as good for predictions made 
using the baseline DES mesh and refined DES mesh. 

Figures 10-13 show comparisons of the predicted midspan pressure coefficients for the three-dimensional 
steady RANS simulations with experimental data for an angle of attack range of 0° to 6° in 2° increments. It 
should be noted that, for each case considered, there were only minimal spanwise variations in the pressure 
distribution. At the lower angles of attack, neither the upper or lower surface pressure distributions agree 
well with experimental data in the first 30% of the chord. As the angle of attack increases, the agreement 
between the lower surface predictions and experimental data improves. The discrepancies between the up- 
per surface predictions and the experimental data increase. This can be explained as follows. At the lower 
angles of attack, the upper and lower surface horns present nearly “equal” disturbances to the flow field. As 
the angle of attack increases, the lower surface horn presents less of a disturbance while the upper surface 
horn presents more of a disturbance. Thus, an increasing angle of attack produces a decreasing region of 
separated flow on the lower surface and an increasing region of separated flow on the upper surface. One 
characteristic that is present in each case is the pressure overshoot that occurs just downstream of the upper 
and lower surface horns. The pressure near the stagnation region appeal's reasonably well-predicted. On 
both the upper and lower surfaces, the flow accelerates around the horn tip, as indicated by the pressure 
drop. A sudden recompression occurs to a pressure level that is too high. The cause of this behavior is not 
currently known. It should be noted that the two-dimensional steady RANS results report by Chi, et al. [21] 
show a different behavior. In the results shown in Ref. [21], better agreement is obtained near the horn. 
The agreement then deteriorates as you move downstream. These results suggest that it is possible that the 
discrepancies shown in Figures 10 through 13 may be due to an overall lack of resolution of the horn tip. 
However, mesh refinement studies by Chung, et al. [32] suggest that the streamwise mesh spacing has only 
minimal impact on solution accuracy for iced wing flow fields. This point is currently being investigated. 
We can also note that, in general, mesh refinement produces some improvement from the coarse mesh to 
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the baseline mesh, when compared to the experimental data. No significant improvement is noted, however, 
after refinement from the baseline mesh to the refined mesh. 

Figures 14-16 show comparisons of the predicted x-component of velocity for the three-dimensional 
steady RANS simulations with experimental data for angles of attack of 0°, 4°, and 6°. Using these images, 
the reattachment location of the primary upper surface flow separation may be estimated by locating the 
position on the chord at which the u=0 contour intersects the upper surface of the wing. At an angle of 
attack of 0°, the reattachment position is well predicted on the coarse, baseline, and refined meshes. As the 
angle of attack is increased to 4°, the downstream extent of the separated region is slightly overpredicted, 
approximately 32% of the chord as compared to 28% of the chord. As the angle of attack is further increased 
to 6°, the extent of the separated region is significantly overpredicted and the upper surface reattachment 
position is shifted downstream until the flow is completely separated. The experimental data shows a flow 
reattachment at approximately 52% of the chord. Although mesh refinement does improve the prediction 
relative to the coarse mesh for both the 4° and 6° angle of attack solutions, no significant improvements 
arc observed when comparing results computed on the baseline and refined meshes. Additionally, at both 
4° and 6°, the transverse extent of the reversed flow region appeal's to be underpredicted (bubble appeal's 
to be thicker in the experimental velocity contours). The extra dissipation present in the coarse mesh is 
evident from the increased spreading of the shear layer in comparison to results predicted using the baseline 
and refined mesh. However, this spreading is not as severe as expected given the relative coarseness of 
the coarse mesh. The seemingly anomalous velocity contours that appeal' in the experimental data just 
downstream of the horn are artifacts from the process employed to generate the contour plots. 

Figures 17-20 show predictions of separation (red) and attachment (blue) locations for an angle of at- 
tack range of 0° to 6° in 2° increments for the coarse RANS mesh, baseline DES mesh, and refined DES 
mesh, respectively. The technique of Kenwright [33] is employed to locate surface mesh elements in which 
a potential separation or attachment line crosses. No attempt is made to reconstruct the actual curves. Ken- 
wright’s technique is based on a phase plane analysis of critical points in the velocity field “close” to the 
surface and, as such, is subject to anomalous results if there is noise in the data. The baseline and refined 
mesh results show fairly well resolved flow features including secondary chordwise separations and corre- 
sponding attachments clearly indicated just down stream of the horn (the closely spaced, roughly parallel 
blue and red curves located just aft of the leading edge). The refined mesh results show a relatively uniform 
primary reattachment across the span in each case. The coarse mesh does not resolve this secondary separa- 
tion at any angle of attack. Additionally, the location of the primary reattachment predicted using the coarse 
mesh does not exhibit the same uniformity that is present in the baseline and refined mesh computations 
except for the 6° case which is fully separated. It is assumed that the coarser, nonuniform surface mesh 
employed in the coarse mesh, as shown in Figure 3(a), is responsible for the appearance of the coarse mesh 
flow fields. Conclusions drawn regarding the attachment locations are consistent with those drawn from the 
velocity contours for the 0°, 4°, and 6° cases. 

To summarize, the steady-state RANS results show an increasing region of separation on the upper 
surface as angle of attack is increased. In all cases, pressure overshoots occur downstream of the upper 
and lower surface horns resulting in an elevated pressure downstream of the horns. In the 4° angle of attack 
results, the downstream extent of the upper surface separation appeal's slightly overpredicted. In the 6° angle 
of attack results, the flow is fully separated. The velocity contour data further suggest that the transverse 
thickness of the separated region is underpredicted. The general trend is that mesh refinement does not 
improve the accuracy of the solution beyond the initial refinement from the coarse mesh to the baseline 
mesh. No significant improvement is observed after further refinement. 


NAS A/CR— 2004-2 13379 


17 




(a) Extruded wing (baseline DES mesh) - lift history 



(b) Extruded wing (baseline DES mesh) - axial force history 


Figure 8: Convergence of RANS solution for extruded wing (baseline DES mesh at a 6° angle of 
attack) 
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RANS SIMULATIONS 



(a) Lift coefficients 


RANS SIMULATIONS 



Angle of Attack (deg) 


(b) Drag coefficients 


Figure 9: Comparison of predicted (RANS) force coefficients and experimentai data for the ex- 
truded wing (coarse RANS mesh, baseline DES mesh, and refined DES mesh) 
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(a) Coarse RANS mesh 


0 deg 



(b) Baseline DES mesh 


0 deg 



(c) Refined DES mesh 


Figure 10: Comparison of predicted (RANS) midspan pressure coefficients and experimental data 
for the extruded wing at 0° angle of attack 
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(a) Coarse RANS mesh 


2 deg 



(b) Baseline DES mesh 


2 deg 



(c) Refined DES mesh 


Figure 11: Comparison of predicted (RANS) midspan pressure coefficients and experimental data 
for the extruded wing at 2° angle of attack 
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(a) Coarse RANS mesh 


4 deg 



(b) Baseline DES mesh 


4 deg 



(c) Refined DES mesh 


Figure 12: Comparison of predicted (RANS) midspan pressure coefficients and experimental data 
for the extruded wing at 4° angle of attack 
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(a) Coarse RANS mesh 
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(b) Baseline DES mesh 


6 deg 
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(c) Refined DES mesh 


Figure 13: Comparison of predicted (RANS) midspan pressure coefficients and experimental data 
for the extruded wing at 6° angle of attack 
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(a) Extruded wing (coarse RANS mesh) 



(b) Extruded wing (baseline DES mesh) 



(c) Extruded wing (refined DES mesh) 


u/U v 



(d) Experimental hot- split- film data 


Figure 14: Comparison of predicted (RANS) midspan u-velocity contours and experimental data 
for the extruded wing at 0° angle of attack 
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(a) Extruded wing (coarse RANS mesh) 



(b) Extruded wing (baseline DES mesh) 



(c) Extruded wing (refined DES mesh) 


uIU" 



xk 


(d) Experimental hot- split- film data [18] 


Figure 15: Comparison of predicted (RANS) midspan u-velocity contours and experimental data 
for the extruded wing at 4° angle of attack 
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(a) Extruded wing (coarse RANS mesh) 



(b) Extruded wing (baseline DES mesh) 



(c) Extruded wing (refined DES mesh) 


u/U v 



(d) Experimental Hot-Split-Film Data [18] 


Figure 16: Comparison of predicted (RANS) midspan u-velocity contours and experimental data 
for the extruded wing at 6° angle of attack 
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(c) Refined DES mesh 


Figure 17: Comparison of predicted (RANS) upper surface separation (red) and attachment (blue) 
locations for the extruded wing at 0° angle of attack 
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(a) Coarse RANS mesh 




(c) Refined DES mesh 


Figure 18: Comparison of predicted (RANS) upper surface separation (red) and attachment (blue) 
locations for the extruded wing at 2° angle of attack 
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(a) Coarse RANS mesh 



(b) Baseline DES mesh 
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(c) Refined DES mesh 


Figure 19: Comparison of predicted (RANS) upper surface separation (red) and attachment (blue) 
locations for the extruded wing at 4° angle of attack 
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(a) Coarse RANS mesh 



(b) Baseline DES mesh 



(c) Refined DES mesh 


Figure 20: Comparison of predicted (RANS) upper surface separation (red) and attachment (blue) 
locations for the extruded wing at 6° angle of attack 
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4.3 DES results 

In this section, results from the two-dimensional and three-dimensional DES computations and com- 
parisons with experimental data [17, 18] arc presented. Results obtained from mesh refinement and time 
step studies arc also included. In addition to integrated quantities such as lift and drag, we present detailed 
comparisons for time-averaged and fluctuating quantities. All DES solutions were initiated from converged 
steady RANS solutions. 

4.3.1 Two-dimensional DES results 

DES computations were performed on the two-dimensional baseline mesh using several different time steps. 
The timesteps chosen were Z\£=0.6ms, Z\f=0.3ms, Z\f=0.15ms, Z\£=0.075ms, and Z\f=0.0375ms at 6° angle 
of attack. DES computations were also performed on the refined mesh at 4° and 6° angles of attack with 
A 1=0.075 ms. Here, At=0. 1 5ms is the time step estimated using the approach described in reference [7] and 
outlined in section 3.5. 

Table 2 shows the time-averaged lift and drag coefficients computed using the baseline DES mesh and 
several different time steps. Note that the lift and drag coefficients appeal - to converge for a time step of 
roughly one-half of the estimated At. Figures 21 and 22 show comparisons between two-dimensional DES 
results and experimental data for the time-averaged x-component of the velocity for angles of attack of 4° 
and 6° for the baseline and the refined mesh. Figures 21 and 22 show that the two-dimensional DES results 
overpredict the size of the separated region downstream of the horn. Figures 21 and 22 and Table 2 indicate 
that the two-dimensional DES results overpredict the force coefficients. Additionally, mesh refinement does 
not have much effect on the the solution. 

Figure 23 shows the instantaneous vorticity plots for A =0.3, A=33.7 and A=67.4. The unsteadiness in 
the flow and the vortex shedding downstream of the airfoil can readily be observed. Figure 23 shows the 
ability of DES method to capture the large-scale vortices for this kind of unsteady flow problems. 

We attribute the significant differences between the simulation results and the experimental data to the 
fact that the simulation is two dimensional. As will be seen in the next section, even though there is no 
spanwise variation in the wing section, the resulting flow field is highly three dimensional. Therefore, we 
would not expect the two-dimensional simulation to accurately model the conditions in the actual flow field. 
Similarly anomalous behaviors were observed in two-dimensional DES results by Kumar and Loth [12], 

Table 2: Two-dimensional DES time-averaged lift and drag coefficients on the baseline mesh, for 
different time steps at 6° angle of attack 


Time Step 

Lift Coefficient (Cl) 

Drag Coefficient ( Cd ) 

Experimental (a = 6.124) 

0.659 

0.105 

Z\f=0.6ms 

0.772 

0.144 

Z\£=0.3ms 

0.869 

0.172 

Z\ £=0.1 5ms 

0.788 

0.154 

At=0. 075ms 

0.877 

0.165 

Z\£=0.0375ms 

0.878 

0.165 
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(a) Two-dimensional baseline mesh 



(b) Two-dimensional refined mesh 







(c) Experimental hot- split- film data [18] 


Figure 21: Comparison of predicted time-averaged u-velocity contours and experimental data for 
the extruded wing (two-dimensional baseline DES mesh, and two-dimensional refined DES mesh 
at 4° angle of attack with At=0. 075ms) 
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(a) Two-dimensional baseline mesh 



(b) Two-dimensional refined mesh 


u/U" 



(c) Experimental hot- split- film data [18] 


Figure 22: Comparison of predicted time-averaged u-velocity contours and experimental data for 
the extruded wing (two-dimensional baseline DES mesh, and two-dimensional refined DES mesh 
at 6° angle of attack with At=0. 075ms) 
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(b) **=33.7 



(c) **=67.4 


Figure 23: Instantaneous vorticity plots for two-dimensional baseline mesh at 6° angle of attack 
with Zlf=0.15ms 
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4.3.2 Three-dimensional DES results: 4° angle of attack, At =0.15ms, baseline mesh 

We now compare DES results obtained using the baseline DES mesh with At- 0. 15ms with the experimental 
data at 4° angle of attack. Table 3 shows a comparison of time-averaged lift and drag coefficients with the 
corresponding RANS results and experimental data. Differences between the steady RANS results and the 
time-averaged DES results are less significant. 

Figure 24 shows comparisons of the predicted time-averaged DES and steady RANS midspan pressure 
distribution with experimental data. The time-averaged DES results show a nearly constant pressure from 
the leading edge to approximately 20% of the chord. This constant pressure region on the upper surface 
indicates that the flow separates from the horn and forms a separation bubble aft of the ice shape. The steady 
RANS predictions show a pronounced pressure decrease in the region ahead of the 10% chord location. 
However, there is no clear constant pressure region. As the flow proceeds aft, the differences between the 
RANS results and the DES results decrease. There are small differences between the DES results and RANS 
results on the lower surface ahead of 30% chord. Aft of this point, the predicted pressures are similar and 
compare favorably to the experimental pressures. 

Figure 25 shows a comparison between the RANS results, the DES results, and the experimental data 
for the time-averaged x-component of the velocity at midspan. The reattachment location on the wing upper 
surface may be estimated by locating the position on the chord at which the u=0 contour intersects the upper 
surface of the wing. From Figure 25, we can estimate that the time-averaged flow is separated for more 
than 40% of chord in the DES simulation and slightly less than 40% of the chord in the RANS simulations. 
The experimental data shows reattachment at approximately the 30% chord position. DES results show a 
secondary recirculation region indicated by the green color ahead of the 10% chord position. This secondary 
recirculation region is much less pronounced in the RANS results. Although the experimental data does not 
show a secondary recirculation region, it is likely that such a region does exist. Regardless, the DES result 
appears to overpredict the extent of this secondary flow. 

Figure 26 shows a comparison of the separation and attachment lines in the predicted flow fields esti- 
mated using the critical point-based approach of Kenwright [33]. The attachment lines shown in Figure 26 
arc consistent with the positions estimated from Figure 25. As noted above, the velocity contours suggest 
a larger secondary recirculation region in the time-averaged DES computations. The broken red spanwise 
line, the separation line for the secondary recirculation, is shifted farther aft in the DES computations. In 
both cases, the blue reattachment line is located very near the horn. There is an anomalous reattachment 
near the trailing edge in the DES results, which is unexplained. 

Figure 27 depicts a comparison between the RMS of the fluctuations in the streamwise velocity compo- 
nent predicted by the DES computation and experimental data. The turbulence intensity was calculated as 
the root-mean-square of the fluctuating x-component of velocity normalized by the freestream velocity. The 
fields arc qualitatively similar, with the exception that the location of the maximum value of the fluctuation 
is shifted upward and the value is larger in magnitude in the predicted results. 

Table 3: Comparison of predicted lift and drag coefficients and experimental data for the extruded 
wing (baseline mesh at 4° angle of attack with At=0. 1 5ms) 



Lift Coefficient (Cl) 

Drag Coefficient (Co) 

Experimental (ct=4.129) 

0.496 

0.0568 

Steady RANS results (baseline mesh) 

0.430 

0.0485 

Time-averaged DES results 
(baseline mesh with Z\t=0.15ms) 

0.414 

0.049 
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Figure 24: Comparison of predicted midspan wing pressure coefficients and experimental data for 
the extruded wing (baseline mesh at 4° angle of attack with Z\ £=0.1 5ms) 
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(a) Steady RANS results (baseline mesh) 



(b) Time-averaged DES results (baseline mesh with A (=0.1 5ms) 
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(c) Experimental hot- split- film data [18] 



Figure 25: Comparison of predicted midspan u-velocity contours and experimental data for the 
extruded wing (baseline mesh at 4° angle of attack with A (=0.1 5ms) 
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(a) Steady RANS results (baseline mesh) 



Figure 26: Comparison of predicted separation (red) and attachment (blue) locations for the ex- 
truded wing (baseline mesh at 4° angle of attack with A (=0.1 5ms) 
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(a) DES results (baseline mesh at Z\t=0.15ms) 
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(b) Experimental hot- split- film data [18] 


Figure 27: Comparison of the RMS of the fluctuations in the u-velocity component for the ex- 
truded wing (baseline mesh at 4° angle of attack with At=0. 1 5ms) 
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4.3.3 Three-dimensional DES results: 6° angle of attack, At =0.15ms, baseline mesh 

We now compare DES and RANS results computed using the baseline DES mesh with At=0 . 1 5 ms with 
experimental data at 6° angle of attack. Table 4 shows a comparison of time-averaged lift and drag coef- 
ficients with corresponding RANS results and experimental data. Lift and drag coefficients from the DES 
computations compare well with the experimental data. However, as we will see below, the predicted flow 
field details do not show the level of agreement suggested by the lift and drag coefficients. 

Figure 28 shows comparisons of predicted steady RANS and time-averaged DES midspan pressures 
with experimental data. The time-averaged DES results again suggest a large separation bubble centered 
just downstream of the midchord position. The predicted bubble position is aft of the position suggested 
by the experimental data. The pressure is seen to be relatively constant from leading edge to approximately 
24% of the chord. Aft of this point, the pressure increases. The DES results show improved agreement 
with the experimental pressures in comparison to the RANS results. There arc no significant differences 
between the DES results and RANS results on the lower surface. The improved agreement between the 
predicted force coefficients and the experimental data shown in Table 4 can be attributed at least in part to 
the strengthening of the recirculation on the wing upper surface relative to the RANS solution. 

Figure 29 shows a comparison between the RANS results, the DES results, and the experimental data for 
the time-averaged x-component of the velocity at midspan. This figure shows that the time-averaged flow 
is separated over nearly the full chord for both the RANS and DES predictions, while the experimental data 
shows a reattachment at approximately 50% chord. As suggested by the pressure data, the time-averaged 
DES velocity data show that the recirculation is shifted aft, as indicated by the region darker blue color 
region. A secondary recirculation region, indicated by the green color adjacent to the wing ahead of the 
15% chord position, is larger in the DES results than in the RANS results, which is consistent with the 4° 
angle of attack case. 

Figure 30 shows a comparison of the separation and attachment lines in the predicted flow fields. The 
flow separates from the horn and reattaches at the trailing edge as indicated by the dark blue line across 
the span. These results arc consistent with Figure 29. Velocity contours from Figure 29 suggest a larger 
secondary recirculation region in the time-averaged DES computations. The broken red spanwise line in 
Figure 30 seems to confirm this result. 

Figure 3 1 depicts a comparison between the RMS of the fluctuations in the x-component of the velocity 
predicted by the DES computations and experimental data. The predicted region of maximum intensity is 
shifted upward and is conspicuously larger in extent than that shown in the experimental data. This suggests 
that the predicted unsteadiness is more pronounced than the physical unsteadiness. However, the region of 
green just downstream of the horn indicates smaller values of the fluctuating component. This implies that 
the secondary recirculation that exists in this region is fairly stable. 

Figure 32 illustrates the unsteady, three-dimensional nature of the DES computation. Figure 32 shows 
isovorticity surfaces (300s~ 1 ) colored by pressure at six different non-dimensional times. The total time 
interval shown in these figures is 0.75s. As can be seen from the images, the flow is highly unsteady and 
there is significant three-dimensionality in the flow once the vortex sheet begins to roll up. The low-pressure 
regions in the pressure field (blue in color) signify the presence of vortex cores. 
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Table 4: Comparison of predicted lift and drag coefficients and experimental data for the extruded 
wing (baseline mesh at 6° angle of attack with At-O. 1 5ms) 



Lift Coefficient (C L ) 

Drag Coefficient ( C D ) 

Experimental ( a = 6.124) 

0.659 

0.105 

Steady RANS results (baseline mesh) 

0.561 

0.083 

Time-averaged DES results 
(baseline mesh with A £=0.1 5ms) 

0.656 

0.098 


6 deg 



Figure 28: Comparison of predicted midspan wing pressure coefficients and experimental data for 
the extruded wing (baseline mesh at 6° angle of attack with Z\ £=0.1 5ms) 
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(a) Steady RANS results (baseline mesh) 



(b) Time-averaged DES results (baseline mesh with At=0. 1 5ms) 
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(c) Experimental hot- split- film data [18] 


Figure 29: Comparison of predicted midspan u-velocity contours and experimental data for the 
extruded wing (baseline mesh at 6° angle of attack with At=0. 1 5ms) 
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(a) Steady RANS results (baseline mesh) 



(b) Time-averaged DES results (baseline mesh with A (=0.1 5ms) 


Figure 30: Comparison of predicted separation (red) and attachment (blue) locations for the ex- 
truded wing (baseline mesh at 6° angle of attack with A (=0.1 5ms) 
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(a) DES results (baseline mesh with Z\i=0.15ms) 
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(b) Experimental hot- split- film data [18] 


Figure 31: Comparison of the RMS of the fluctuations in the u-velocity component for the ex- 
truded wing (baseline mesh at 6° angle of attack with At=0. 1 5ms) 
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(c) r=47.2 


(e) £*=60.7 



(d) t*= 53.9 



(f) £*=67.4 




Figure 32: Isovorticity contours (300s 1 ) colored by pressure for the DES baseline mesh at 6° 
angle of attack with Z\ £=0.1 5 ms 
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4.3.4 Three-dimensional DES results: 6° angle of attack, At =0.075ms, baseline mesh 

In this section, results obtained using the three-dimensional baseline DES mesh with Af=0.075ms are com- 
pared with RANS results and experimental data for a 6° angle of attack. Table 5 shows a comparison of 
time-averaged lift and drag coefficients with the corresponding RANS solutions and the experimental data. 
In general, the lift and drag coefficients compare well with experimental data. 

Figure 33 shows comparisons of predicted time-averaged midspan pressures with experimental data. 
Overall, the time-averaged DES results at A (=0.075 ms arc quite similar - to the A (=0.1 5ms results. As 
before, the agreement between the predicted DES results and the experimental data is not commensurate 
with that suggested by the lift and drag data. Figures 34 and 35 compare the time-averaged x-component of 
the velocity at midspan and the separation and attachment lines in the predicted flow fields. Again, the extent 
of the separated region is overpredicted. The separation and attachment line images are consistent with the 
velocity contours. Figures 34 and 35 again suggest the presence of a significant secondary recirculation 
region. Figure 36 compares the predicted turbulence intensities with experimental data. This image shows 
that the intensity of the fluctuations is increased in the predictions and occurs over a larger region than in 
the experimental data. The stability of the secondary flow structure is again evident. 

Figure 37 shows isovorticity surfaces (300s -1 ) colored by pressure at six different non-dimensional 
times. The total time lapse shown in these figures is, again, 0.75s. The low-pressure regions in the pressure 
field (blue in color) indicate the presence of vortex cores. These images indicate that the flow field is highly 
unsteady and three dimensional. 

Table 5: Comparison of predicted lift and drag coefficients and experimental data for the extruded 
wing (baseline mesh at 6° angle of attack with At =0.07 5 ms) 



Lift Coefficient (Cl) 

Drag Coefficient (Co) 

Experimental (a = 6.124) 

0.659 

0.105 

Steady RANS results (baseline mesh) 

0.561 

0.083 

Time-averaged DES results 
(baseline mesh with A (=0. 075ms) 

0.677 

0.098 
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Figure 33: Comparison of predicted midspan wing pressure coefficients and experimental data for 
the extruded wing (baseline mesh at 6° angle of attack with ZU=0.075ms) 
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(a) Steady RANS results (baseline mesh) 



(b) Time-averaged DES results (baseline mesh with At=0. 075ms) 
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(c) Experimental hot- split- film data [18] 


Figure 34: Comparison of predicted midspan u-velocity contours and experimental data for the 
extruded wing (baseline mesh at 6° angle of attack with At=0. 075ms) 
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(a) Steady RANS results (baseline mesh) 



(b) Time-averaged DES results (baseline mesh with At - 0.075ms) 


Figure 35: Comparison of predicted separation (red) and attachment (blue) locations for the ex- 
truded wing (baseline mesh at 6° angle of attack with A (=0.075 ms) 
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(a) DES results (baseline mesh with Z\t=0.075ms) 
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(b) Experimental hot- split- film data [18] 


Figure 36: Comparison of the RMS of the fluctuations in the u-velocity component for the ex- 
truded wing (baseline mesh at 6° angle of attack with /\f =0.075 ms) 
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(a) r=33.7 


(c) r=47.2 


(e) £*=60.7 


(b) £*=40.4 



(d) £*=53.9 



(f) £*=67.4 





Figure 37: Isovorticity contours (300s 1 ) colored by pressure for the DES baseline mesh at 6° 
angle of attack with At=0. 075ms 
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4.3.5 Three-dimensional DES results: 6° angle of attack, At =0.075ms, refined mesh 

In this section, DES results obtained using the refined DES mesh with Z\t =0.075 ms are compared steady 
RANS results and experimental data for a 6° angle of attack. Table 6 shows a comparison of time-averaged 
lift and drag coefficients with the corresponding RANS results and experimental data. The lift and drag 
coefficients again compare well with the experimental data. 

Figure 38 shows that the constant pressure region terminates around x j c ~ 0.28. Even in this case, 
the extent of separation bubble is overpredicted when compared to experimental data. The lower surface 
pressure distribution is similar for the RANS and DES results. Figure 39 shows comparisons between 
predicted time-averaged x-component of the velocity at midspan and experimental data. From this figure we 
can note that even though the separation region extends up to the trailing edge, the intensity of the reverse 
flow is decreased as the boundary layer reattaches. This is evident in the decrease of the blue region near the 
trailing edge. The secondary recirculation region is again evident. Separation and reattachment lines from 
Figure 40 indicate flow features consistent with those shown in Figure 39. 

The turbulence intensity contours for the extruded wing in Figure 41 agree qualitatively with the corre- 
sponding experimental data. Here the location of maximum intensity is shifted slightly upward and to the 
right. Isovorticity surfaces (300.s 1 ) colored by pressure at six different non-dimensional times are shown in 
Figure 42. Here, the vortex cores arc more readily discerned than in the other results. The improved spatial 
resolution is most likely responsible for this result. 

Table 6: Comparison of predicted lift and drag coefficients and experimental data for the extruded 
wing (refined mesh at 6° angle of attack with At=0. 075ms) 



Lift Coefficient (Cl) 

Drag Coefficient (Co) 

Experimental (a=6. 1 24) 

0.659 

0.105 

Steady RANS results (refined mesh) 

0.558 

0.085 

Time-averaged DES results 
(refined mesh with Z\f=0.075ms) 

0.666 

0.090 
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Figure 38: Comparison of predicted midspan wing pressure coefficients and experimental data for 
the extruded wing (refined mesh at 6° angle of attack with ZU=0.075ms) 
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(a) Steady RANS results (refined mesh) 



(b) Time-averaged DES results (refined mesh with Z\t=0.075ms) 
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(c) Experimental hot- split- film data [18] 


Figure 39: Comparison of predicted midspan u-velocity contours and experimental data for the 
extruded wing (refined mesh at 6° angle of attack with ZU=0. 075ms) 
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(a) Steady RANS results (refined mesh) 



(b) Time-averaged DES results (refined mesh with Z\t=0.075ms) 


Figure 40: Comparison of predicted separation (red) and attachment (blue) locations for the ex- 
truded wing (refined mesh at 6° angle of attack with Z\ /:=(). 075ms) 
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(a) DES results (refined mesh with ZU=0.075ms) 
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(b) Experimental hot- split- film data [18] 


Figure 41: Comparison of the RMS of the fluctuations in the u-velocity component for the ex- 
truded wing (refined mesh at 6° angle of attack with Z\ /:=(). 075ms) 
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(a) t*= 33.7 



(c) t*= 47.2 



(e) £*=60.7 



£*=40.4 



(d) £*=53.9 



(f) £*=67.4 


Figure 42: 

Z\t=0.075ms 


Isovorticity contours (300s ) colored by pressure at 6° angle of attack with 
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4.3.6 Mesh refinement and time step refinement comparisons 

This section su mm arizes the effects of mesh refinement and timestep refinement for three-dimensional DES 
computations. Table 7 and Figures 43(a) and 43(b) show comparisons of force coefficients from steady 
RANS and DES solutions. In general, we can see that there is no significant improvement in the accuracy 
of the solution after mesh refinement. There is some improvement in the lift prediction for the DES solution 
for the 4° angle of attack case with time step refinement. Figure 44 shows comparisons of the predicted 
time-averaged, three-dimensional midspan pressure coefficients for the extruded wing with experimental 
data using the baseline mesh and the refined mesh at an angle of attack of 6°. From the figure and Table 7, 
we can note that mesh refinement and time-step refinement do not significantly improve the accuracy in the 
solution when compared to the experimental data. This suggests that the baseline mesh with At=0.15ms, as 
estimated from the guidelines given by Spalart [7], is sufficient for this case. The above conclusion is further 
supported by the stream wise velocity component images (Figures 29, 34 and 39), separation and attachment 
line images (Figures 30, 35 and 40), and turbulence intensity images (Figures 31, 36 and 41). 

The isovorticity contour images (Figures 32, 37 and 42) confirm that there is significant three-dimension- 
ality in the flow. As the timestep and mesh are refined, we also note that the intensity of the reverse flow 
is decreased as the boundary layer reattaches. This is evident from Figures 34 and 39 through a decrease 
in the extent of the blue region near the trailing edge. The turbulence intensity for the refined mesh with 
A f =0.075 ms at a 6° angle of attack compares better with the corresponding experimental data. This can be 
seen from Figure 4 1 . 

To validate the three-dimensional DES solutions obtained for the baseline and refined meshes, the pe- 
riodogram function of MATFAB [34] was used to analyze the time histories of the wing lift coefficient. 
Figure 45 shows the lift force history for the refined mesh with At =0.075 ms at a 6° angle of attack. Fig- 
ure 45(a) shows the lift history for both RANS and DES computations, whereas Figure 45(b) shows only 
the lift history for the DES computation. It is representative of the other simulations. The power spectral 
density (PSD) of the time history of the lift coefficient is computed and plotted against the frequency of the 
signal. The time interval was 0.75s to 1.5s after the initiation of the DES computation. For At=0.075ms, 
this corresponds to 10,000 time steps. For Af=0.15ms, this corresponds to 5,000 time steps. 

Figure 46 shows a comparison of the PSD signatures for the three-dimensional DES solutions on the 
baseline and refined meshes with Af=0.075ms for 6° angle of attack. From the figure, it is apparent that 
similar signatures arc obtained up to 15Hz. Figure 46(b), which is the zoomed version of Figure 46(a) in 
the frequency range 10-100Hz, shows differences above 15Hz. Figure 47 shows a comparison of three- 
dimensional DES solution on baseline mesh with Af=0.15ms and A f =0.075 ms for 6° angle of attack. 
Again, Figure 47 shows very similar signatures up to 10Hz. Figure 47(b), which is the zoomed version 
of Figure 47(a) in the frequency range of 10-100Hz, again shows there arc some quantitative differences at 
frequencies of approximately 12Hz. 

These results indicate that both the time-step refinement and mesh refinement only significantly affected 
the signal at frequencies higher than 12-15Hz. The power at these frequencies is at least an order of mag- 
nitude less (range 10 -04 to 10 -06 ) than the power at the lower frequencies. Taken together with the fact 
that the time-averaged quantities are so similar, we can conclude that these refinements did not significantly 
affect the results of the simulations. Since the character of the solutions did not change with these refine- 
ments, this suggests that we have valid DES results within the context of the AVUS flow solver. This does 
not address the accuracy of the simulation. 
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Table 7: Comparison of lift and drag coefficients for all three-dimensional RANS and DES results 
(coarse RANS mesh, baseline DES mesh and refined DES mesh) 



C L 

Angle of Attack ( a ) 

0 

2 

4 

6 

Experimental (a) 

-0.024 

2.1 

4.129 

6.124 

Experimental 

0.063 

0.289 

0.496 

0.660 

3D-RANS (Coarse) 

0.065 

0.265 

0.443 

0.556 

3D-RANS (Baseline) 

0.065 

0.251 

0.430 

0.561 

3D-RANS (Refined) 

0.066 

0.254 

0.438 

0.558 

3D-DES (Baseline, Z\/=0.15ms) 

N/A 

N/A 

0.414 

0.656 

3D-DES (Baseline, Z\/=0.075ms) 

N/A 

N/A 

0.438 

0.670 

3D-DES (Refined, Z\/=0.075ms) 

N/A 

N/A 

0.445 

0.660 



c D 

Angle of Attack (a) 

0 

2 

4 

6 

Experimental (cr) 

-0.024 

2.1 

4.129 

6.124 

Experimental 

0.042 

0.041 

0.076 

0.105 

3D-RANS (Coarse) 

0.039 

0.042 

0.059 

0.098 

3D-RANS (Baseline) 

0.036 

0.037 

0.048 

0.083 

3D-RANS (Refined) 

0.035 

0.035 

0.048 

0.085 

3D-DES (Baseline, Z\/=0.15ms) 

N/A 

N/A 

0.049 

0.098 

3D-DES (Baseline, Z\f=0.075ms) 

N/A 

N/A 

0.050 

0.098 

3D-DES (Refined, Z\f=0.075ms) 

N/A 

N/A 

0.047 

0.090 
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THREE-DIMENSIONAL RANS AND DES SIMULATIONS 



Angle of Attack (deg) 


(a) Lift coefficients 


THREE-DIMENSIONAL RANS AND DES SIMULATIONS 



Angle of Attack (deg) 


(b) Drag coefficients 


Figure 43: Comparison of predicted wing lift and drag coefficients and experimental data for the 
extruded wing (coarse RANS mesh, baseline DES mesh and refined DES mesh) 


NASA/CR— 2004-213379 


60 


6 deg 


-3 



-0.1 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 


x/c 


Figure 44: Comparison of predicted midspan wing pressure coefficients and experimental data for 
the extruded wing (coarse RANS mesh, baseline DES mesh and refined DES mesh at 6° angle of 
attack) 
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(a) Lift history for RANS and DES computations 

3D-DES (refined mesh, dt = 0.075ms) 
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(b) Zoomed version of Figure 45(a) showing DES lift history 
Figure 45: Lift history for refined mesh with At=0. 075ms 



3D-DES (refined mesh, dt = 0.075ms) 
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(a) Power spectral density 



(b) Power spectral density (zoomed version of Figure 46(a) 


Figure 46: Comparison of power spectral density plots for baseline mesh and refined mesh at 6° 
angle of attack, with Z\t=0.075ms 
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(b) Power spectral density (zoomed version of Figure 47(a) 


Figure 47: Comparison of power spectral density plots for baseline mesh at 6° angle of attack 
with, ZH=0.15ms and Zlt=0.075ms 
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5. Assessment of Existing Capabilities and Recommendations 

In this section, we provide an assessment of the technologies employed in this effort and make a few recom- 
mendations regarding potential roles they may play in icing effects analysis. 

5. 1 Tool assessmen t 

1. The mesh generation tools employed here performed adequately for the iced wing problem. Solid- 
Mesh proved somewhat difficult to use due to the fact that control of point spacing in the interior of 
the domain was challenging. The source-based approach employed in the GridTool/VGridns combi- 
nation provided an effective mechanism to control the point spacing in the domain. However, we did 
encounter some difficulties generating a complete mesh. These problems were solved by trial-and- 
error adjustment of the source strengths. 

2. Effective utilization of the AVUS flow solver was challenging because we were unable to consistently 
execute the program on more than 64 processors. Considering the time required to obtain DES so- 
lutions, this represents a significant impediment to routine deployment of the flow solver for DES 
computations - even if sufficient resources were available. Additionally, while analysis of the re- 
sults indicates that the DES solutions arc valid, questions remain about the quality of the solution. 
In particular, if additional resolution is required near the horn, it may make the computation cost 
prohibitive. 

5.2 Recommendations 

1. For constant section unswept wings, results obtained here suggest that there is no need to employ 
three-dimensional steady RANS simulations. Although there is some three-dimensional flow present 
in the solutions, it seems to have little influence. Additionally, the steady three-dimensional RANS 
simulations do not model the three-dimensionality present in the unsteady DES computations. Two- 
dimensional steady RANS computations appear to be adequate for configurations of this type. 

2. Conversely, there is no compelling reason to employ two-dimensional DES computations. As re- 
ported here and elsewhere [12], the two-dimensional DES computations significantly overpredict the 
effects of the unsteady flow field. This is an inherent limitation of the two-dimensional approxima- 
tion. 

3. It is our opinion that DES should continue to be investigated for application to the icing effects 
problem. Other cases should be simulated so that a better assessment of DES capabilities may be 
obtained. However, the results observed here indicate that only configurations for which detailed 
experimental flow field data is available should be considered. 

4. It is also our opinion that, due to the significant unsteadiness observed in the experimental results [35], 
the extruded GLC305/944-ice shape may represent a particularly difficult case. Therefore, computa- 
tions should be performed for configurations with less challenging horn ice accretions. Because of 
the detailed LDV measurements, one candidate configuration is the swept wing with simulated ice 
studied by Bragg, et al. [36, 37]. 

5. Because we believe that the hybrid RANS/LES strategy is the best near-term option for computing 
these complex flow fields, other hybrid RANS/LES implementations should be investigated along 
with other flow solvers to address questions about AVUS-specific behaviors. 
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6. Summary 

The problem of flow field simulation for iced wing configurations is a complex one that severely taxes ex- 
isting capabilities for geometry modeling, mesh generation, and flow solution. In this report, we focused on 
activities associated with the first year of an effort to investigate the capability of existing flow simulation 
methodologies to predict the separated, highly unsteady flow fields associated with wings with significant 
horn-ice accretions. The effectiveness of Detached-Eddy Simulation (DES) as a tool for predicting icing 
effects was evaluated. The AVUS code was employed to compute solutions for an iced wing configuration 
using DES and steady Reynolds Averaged Navier-Stokes (RANS) equation methodologies. The configura- 
tion was an extruded GLC305/944-ice shape section with a rectangular planform. The model was mounted 
between two walls so no tip effects were considered. The numerical results were validated by comparison 
with experimental data for the same configuration. 

We presented the results of three-dimensional steady RANS computations and two- and three-dimension- 
al DES computations performed for the extruded wing. The time-averaged DES computations showed im- 
provement in lift and drag results near stall when compared to steady RANS results. However, flow field 
details did not show the level of agreement suggested by the integrated quantities. 

The benefits of employing three-dimensional DES computations are not clear at this time. While com- 
parisons of integrated time-averaged quantities such as lift and drag with experimental data arc improved 
relative to steady RANS predictions, commensurate agreement was not obtained for detailed flow field quan- 
tities. In particular, DES results showed more extensive flow separation than the experimental data. Further, 
the DES results showed a significant secondary recirculation that was not present in the experimental data. 

Based on our results, we believe that DES may prove useful in a limited sense to provide analysis of iced 
wing configurations when there is significant flow separation, e.g., near stall, where steady RANS computa- 
tions arc demonstrably ineffective. Additionally, we believe that DES may prove useful for flow fields with 
extensive three-dimensionality such as swept finite wings or wings with genuinely three-dimensional ice 
accretions. However, more validation is needed to determine what role DES can play as paid of an overall 
icing effects prediction strategy. 
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